clc;
warning('off');
M1=5;
K=2;
Ne=2;
R=1;
alpha=0.03;
e=zeros(1,K);
x_axis=5:1:5;
y_axis1=zeros(1,length(x_axis));
y_axis2=zeros(1,length(x_axis));
% y_axis_SDR=zeros(1,length(x_axis));
num=100;
a=0;
for i=1:num
    sign=1;
    y_axis_temp1=zeros(1,length(x_axis));
    y_axis_temp2=zeros(1,length(x_axis));
    disp(i);
    for Nr=x_axis
        hr=(randn(1,Nr)+1i*randn(1,Nr))/sqrt(2);
        G=(randn(Nr,M)+1i*randn(Nr,M))/sqrt(2);
        Dr=(randn(Ne,Nr,K)+1i*randn(Ne,Nr,K))/sqrt(2);
        disp(Nr);
        for k=1:K
            e(1,k)=alpha*sqrt(norm(Dr(:,:,k),2)^2);
        end
       try
        [ w,F,power1 ] = Alternate_optimization( hr,G,Dr,e,R,M,Nr,Ne,K );
        [ w_SDR,F_SDR,power_SDR ] = SDR( hr,G,Dr,e,R,M,Nr,Ne,K );
        catch
            a=a+1;
            sign=0;
            disp('false')
            break;
        end
        %y_axis_random_phase_temp(1,int8(alpha/0.04)+1)=power_random;
        y_axis_temp1(1,int8(Nr/5))=power1;
        y_axis_temp2(1,int8(Nr/5))=power_SDR;
    end
   if sign
       y_axis1=y_axis1+y_axis_temp1;
       y_axis2=y_axis2+y_axis_temp2;
   end
end
y_axis1=y_axis1/(num-a)
y_axis2=y_axis2/(num-a)
y_axis1=10*log10(y_axis1*1000)
y_axis2=10*log10(y_axis2*1000)
%plot(x_axis,y_axis,'c->','LineWidth',1.5,'MarkerSize',6);
% hold on;
% plot(x_axis,y_axis_SDR,'r-*','LineWidth',1.5,'MarkerSize',6);
%hold on;
plot(x_axis,y_axis1,'b-o','LineWidth',1.5,'MarkerSize',6);
hold on;
plot(x_axis,y_axis2,'r->','LineWidth',1.5,'MarkerSize',6);
hold on;
grid on;
xlabel('alpha'); 
ylabel('Worst-case Transmission power (dbm)');
legend('Alternate_optimization M=7','Alternate_optimization M=9','Fix Phase M=7','Fix Phase M=9')